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Gerald E. Wilson 

Laboratory for Theoretical Studies 
ABSTRACT 

A comparison is made among three methods of obtaining highly accurate 
approximate expressions for the average motion of plasmas in slowly varying 
electric and magnetic fields. It is shown that no inconsistencies are revealed 
by this comparison and that, in fact, two of the techniques are completely equi- 
valent while most of the results of the third can be obtained from the equiva- 
lent pair. As a byproduct of the comparison, a function is obtained which can 
be shown to be conserved along the guiding center trajectory. 
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A COMPARISON OF THREE SETS OF HIGHER ORDER 
ADIABATIC PLASMA EQUATIONS 

I. INTRODUCTION 

Three papers Frieman et al . 4 (I), Macmahon^ (II), and Wilson^ (HI) have 
been presented describing separate approaches to the task of obtaining highly 
accurate approximate solutions to the Boltzmann and/or hydromagnetics equa- 
tions for slowly varying fields of arbitrary geometry. Techniques (I) and (III) 
are very closely related as a closed system of self-consistent equations, cor- 
rect to first order in € = m/e the particle mass to charge ratio, is obtained by 
maintaining both microscopic and macroscopic equations. In (II), the e expansion 
is introduced directly into the heirarchy of moments equations obtained from 
the Boltzmann equation. First-order portions of certain quantities can be de- 
termined explicitly in terms of higher moments but for this totally macro- 
scopic description, closure may be obtained only through the imposition of 
additional approximations . 4 

It may easily be verified by inspection of the original papers that the assump- 
tions upon which these works are based are completely equivalent in most plasma 
density and temperature realms of interest. Each set of conditions may be re- 
duced to assuming that the expansion parameter e may be as small as desired. 
Moreover, (I) and (II) are based on the Boltzmann equation with the Lorentz 
force expression while (IQ) is based on a separate solution of the Lorentz 
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equation combined with an expression of phase space conservation; an equivalent 
approach. Thus in as much as comparable expressions are calculated the re- 
sults of the three must reasonably be the same. However, the veracity of this 
supposed equivalence is far from obvious, even though one of the results of (II) 
is verified in (HI) . 5 The equivalence of the microscopic equations in (I) and (HI) 
is totally obscured by the nature of the special variables chosen in (HI), and 
other results of (H) which should be obtainable from (HI) have not been calcu- 
lated. The purpose oi this paper is to show that the expected equivalence does 
obtain; that the additional results of (H) can be derived from (HI), and that (I) 
and (HI) are as alike as to allow the calculation of the results of one from those 
of the other. A by-product of the inspection is the discovery of the form of a 
function which is conserved along the higher order guiding center trajectory. 

With this in mind, in Section H the relevant results of (I) are rederived in 
the notation of (HI), while Section HI supplies the analogous results of (HI). 

In Section IV, a detailed comparison is made between (I) and (HI), the results of 
(H) are verified by use of (HI), and the conservation of a proposed function along 
the guiding center path is demonstrated. 

H. HIGHER ORDER CHEW- GOLDBERGER- LOW THEORY 

(I) is an extension to higher order in e of a procedure proposed by Chew, 
Goldberger and Low . 6 The method is based on an ordering scheme that re- 
quires, in addition to the assumption of slowly varying fields, plasma densities 
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and temperatures such that the plasma frequency is of order of magnitude e* 1 
and the Debye length of order e. Implied in this (as in all three theories dis- 
cussed in this paper) is the assumption that as e - 0 the component of the electric 
field parallel to B also vanishes. 

One proceeds, then, to solve the collisionless Boltzmann equation 

+ (v V) f + 7(E 1 +E || L + vxB) • V y f = 0 (1) 

in the limit as e.- 0. In (1), the subscripts ii and i refer to the directions parallel 
to and perpendicular to the magnetic field B; V y is the divergence operator in 
velocity space; L = B/| b| with M in the direction of the principle radius of 
curvature of the lines of magnetic induction and N = L x M . 

Referring to Figure 1, to all orders, v may be represented by 

v = U + rj£ + a(M sin v + N cos v) (2) 


where 


(a) U = 


ExB 

B 2 


(b) ri - L ■ v 


(C) CT - Iv^-Uj 3 ( v ' | 


(d) V 



(i) 
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Figure 1. Cylindrical Coordinate Sytlem (Different From That of (I)) 

The distribution function f is assumed to be represented by a power series 
in e and each term f i is given in general by f t = 7. + where 



Under the conditions stated, it follows from (1) that the zero order term of 
the distribution Auction is independent of v. This knowledge plus the assump- 
tion that f j is periodic in v makes it possible to arrive at an equation for 7 0 by 
averaging (1) over v and neglecting all but the lowest order terms. 

To next order, f is not independent of v and the solution is algebraically 
more difficult to obtain. One may invert (1) to find dl l /dv in terms of the 
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previously defined f 0 . By integrating df x fdv over phase from zero to v and 
imposing 4b, one may obtain f v The result is: 



where: 

I. s NN “ MM , 
1^ s Mfi + NM , 

V = U + vL . 


The notation d/dt refers to the convective derivative d/dt + (U + T 7 L) • V and 
the double dot tensor product notation is exemplified by ab .* V* = • • (b- V) c. 

By substituting (5) in (1) and integrating over v one may obtain an equa- 
tion of the following form for the first two terms of the power series 



I 


1 

1 



for f: 


where 


|f + (“ * V) f + 0 


S) 


T# ■ • 


(6a) 


e /v dV 6 o' 1 2 ~ cx x ~ 

a = V+fLx + "g ~2B (LxV)B + e 2B ^ : 


(6b) 


p = - -y i A : VV 


+ e 


2B 


RV . f u dV dvl 

B \ M B ' dt " N B * dfc J 


_ dV _ ~~ dV 

+ NL : VLN • ^ - N L : VLM • g*- 


1 

2(n* £uM : VL - M • gf UN : Vl) 


(6c) 
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where: 

I y = NM - MN , 

I A = MM + NN . 

To obtain a closed set of equations suitable for the study of plasma behavior, 
one must add to (6) equations which will determine the fields E and B self- 
consistently. This has been accomplished and reported in the original paper 8 
and will not be discussed here. Equations (5) and (6) are the relevant equations 
to be studied in Section IV. It will be observed in Section m that an equation 
which serves the same purpose as (6) but which has a radically different ap- 
pearance can be obtained. We will show that the two equations are consistent. 
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in. THE GUIDING CENTER PLASMA 

The approach first outlined by Grad** » 10 and extended in (HI) is quite dif- 
ferent from that of the previous section. It is based on guiding center orbits 
obtained by solving the single particle equation of motion: 

3 v 1 

ST = T [E,+e,£+.xB] . ( 7 ) 

As outlined in (HI), 11 the higher order terms of this solution can be calculated 
in terms of a "nice" 1 ^ se t of variables chosen to have the property that their 
equations of motion are free of any rapidly oscillating terms depending on phase. 
If the particle coordinates are represented in terms of "nice" variables as (P, 

2 , H, $>) these variables are related to those of Section n by the following 

P = r-gLxv' (8a) 


2 2 

2 



[ctM+Lxv'] • 


dV 

dt 


e a * r 

+ -g ~2 I a sin v cos v + Ip s 


in 2 : VV (8b) 


H 


dL 


- 7j - -g- U [ctM + Lxt'] : V-L - g [er M + Lx v'] • ■jjj: 


e ~ 2 q [l a sin v cos v + I * sin 2 v] : VL (8c) 


8 


V + 0(e) 


(8d) 


<t> - 


where the order e(0(e)) terms of (8d) are known bit never used. 

In the present case, the variables may reasonably be referred to as the 
guiding center variables since P is the position of the guiding center and the 
motion of the guiding center may be completely described in terms of P, Z, H 
and their derivatives. One may also describe the instantaneous motion of the 
particle itself in terms of these variables by using the inversion relations 

v = U(r) + 7jL(r) + cr {m (r) sin v + N(r) cos v} (9a) 

r = P + g7~y (N(r) sin v -jfi(r) cos v) . (®b) 


In (9), the lengthy transformations (8) must be used to develop the expressions 
for r and v as power series expansions to the order of accuracy desired. ^ 
The time derivatives of (8) are: 


SP 

St 


U + HL + -§ L 


x 


dV 
dt + 


el 2 ~ 

B 2B (k x B 


B 1 


^ + UM: Vl) + ^ (l y + -§ l fi ) : VL (10a) 
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= - ~2 I x : V (U + H L) 


s 2 

_ JM ~ 

dV 

N ~ 

dV] 

+ e 2B 

BV -jg N * 

dt 

B M * 

dt j 

A A 

~ * dV 

A A 

Aft 

dV 

+ ML : 

VLN • jt 

- NL 

: VLM 

' dt 


T (« 


S*f : ’C 


/- dL 
( N - dt 


UM : VL- M • ^ UN : VLJ 


+ T^ : VVI i : 7V 'Idl ^ : VV 

- • 


/ E ll c£\ - dL 3 - dV „„ d . dV 

2 \ e + U ' dt ) M ' dt + 2 M ’ dt l x : V V + dt M ' dt 


-y 3 u 

+ -fir [ML : VU- UM : V-£] I A : V-L+ € ML : VL I x : VL (10b) 


E il dL l 2 „ 

T +u 'dF + T v ' L 


IT BV-{|(»T*^+UM:Vl)-|(n. ^+UN:Vl)} 
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+ MM : VUIj : VL - MM : V-LI^ : VU 


_ 1 


1 d 




M- 

W:: 


ff:- 

M&- 

1 

t v -, 


& 






2 V * VI* :VL- jdt ** : VL 


el 

B 


jlJL: Vl(m • 37 + UM: Vl) 


1 dL -\ „ / E || d£\ 

“ 2 l i : VlJ V M * dt + UM : VL ) ‘ ML : VL \T + U * gt j 


_d_ 

dt 


/ ^ dL „ dU -~| 

|M • dF + U M : VL j - M • dT V ‘ L J 


e £ 3 

V • L - g -y M L : VLV • L 


HI (3 - dL 1 ~ 

- e g [2 M * dt + 2 ^ ** : 7 L J 

e L dl) / dL „ dU /- dL . Si 

+ bK dt (« • d F+ UN: VL] - N • ^ (m • jT + Vl)J 

v> 

+ X • 3F UN : VL- N • UM : VlJ (10c) 


8$ B dM 1 , 

2v Jt = T +N 'dT + T I y :VV - 

Here d/dt and V are defined as d/d t + (U + HL) • V and U + HL respectively. 

These definitions agree with those of the previous section to lowest order. In 

the next section the notation V and d/dt will be used only when higher order terms 

are to be neglected. 
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In (HI), a further change of variables is constructed so that one velocity 
space component is characterized by the magnetic moment series correct to 
order e. The resulting equation is of a simpler form and is more useful for 
the purposes of (in). This equation will be discussed and compared with an 
analogous equation in Section IV, but at this point we choose to present the eon- 
sequences of the results thus far obtained. That is, one may determine the 
probability density F of guiding centers in the volume element d 3 PldIdHd$ 
through the conservation relation 



The theorem of phase independence of Reference 12^ and Equation (10d) jus- 
tify the absence of a term involving 8<l>/8t in (11). 

Equation (11) may be written in a form more suitable for comparison with 
(6) as follows 
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Since the right-hand- side of (12) does not vanish, this equation expresses the 
fact that F is not conserved along a guiding center path. We will note presently 
that a function can be found that does satisfy such a criterion. Again, the dis- 
cussion of the self-consistent closure of the system is to be omitted since (8), 
(10), and (12) are the results necessary for comparison with (6). 


IV. COMPARISON 

With the notation established and all necessary results stated, we can now 
proceed to demonstrate that the results of Sections II and HI are consistent with 
each other. The connection with (II) will be dealt with separately since the 
microscopic equations of (I) and (ID) are of primary interest. 

The two distribution functions of Sections II and III 


I — ^ 

\t, 2 . T7. vj 


and 


(p, y , H, 4>) 


may be related by equating the number of particles in d 3 t<t d<r d t? d v to the 
number of guiding centers in d 3 P 2 d 2d H d an equation expected to be true 
on physical grounds. That is 


f d 3 r<rd<xdi?d v = Fd 3 P2d2dHd«. 


(13) 
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The transformation Equations (8) make It possible to relate the particle 
and guiding center volume elements through the Jacobian of this transformation 
so that 


d 3 radad7]dv - Jd 3 P£dXdHd4> . 


(14) 


From (13) and (14) one may conclude that 



(15) 


In (15) the argument $ has been dropped in F because, as previously pointed 
out, F is known to be phase independent. 

In an attempt to arrange (15) so that both sides of the equation are in terms 
of the same variables it is advantageous to expand f in Taylor series about the 
guiding center variables as follows 



fp, \ , H. *J + (r-P) • Vf 

(a 2 1 2 \ df df 

+ [- 2-^2 j ^\ +(77 " H) ^ + (V_ 


df 

®) T 7 ‘ ( 16 ) 


Recalling from Section n that one may always write f as 


f 



l 2 

2 


, H, 



2 



(17) 
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we can substitute (5) into (17) and then (17) and (3) into (16) to obtain 



It is important to understand exactly what has been done here. Equation 
(17) represents no transformation of variables but only a renaming of arguments. 
The expression for ? is obtained by merely insertijg into (S) the arguments in- 
dicated in (17). The function T(r, >-- 3 /2, 17) is determined by (6) with a, /3, and r 
written as functions of <r 2 /2 and v- In like manner, the function f (P, £ 3 /2,H) 
that appears on the right-hand-side of (18) is determined by (6), provided the 
arguments of a, ft and T are now replaced by P, X 3 /2 and H. From this point on 
then we assume (18) to have been enforced and all expressions of both Sections 
II and m are to be considered functions of P, X 3 /2 and H. 

Taking note that the right hand side of (18) is independent of and recall- 
ing that F is also independent of $ one may conclude that J must be phase in- 
dependent for (15) to hold. From (8) it can be shown that this is true; J is 


independent of <f and of the form 


1 - 


£ r> 

bT m 


dV <■* /\ 

dF + "b' ml 


VL - 


is h 


vv 


(19) 


thus revealing no inconsistency. 

Proceeding, one may invert (15) to obtain 

T = J" 1 F - e y . (20) 


Insertion of (20) and (6) yields 


d' F 
dt 


Fd'J 

J dt 


t d ' y 

e J ~ST 


( 21 ) 


where: 



Using some shorthand notation we rewrite (12) as 


DF 

W 


+ FA 


0 


( 22 ) 


where: 


DF _ dF ( Sk* „\ _ S \ 2 / <JF 8H dF 

W = Tt r \jS t * V / F + St "M\ + St aff 

d W 
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and 


A 


V 



A |H 
<?H 8t 


Subtracting (21) from (22) one may obtain the following relation which must 
be satisfied by Equations (6) and (10) 




r ) 


dF 

W 


F fd' J 

= “ y ~3t + J A “ 
_ _ 


ej 


d'y 

~St' 


(23) 


Equations (10), (19), and (21) yield 


* = 



(24) 


which is to be expected since the left-hand-side of (23) contains no term propor- 
tional to F. 

The comparison expressed by (23) now acquires additional significance . 

One may note that y is obtained by an entirely separate calculation from (6) and 
(10), and that d' y/dt need be computed to lowest order only. Thus a simple cal- 
culation yields a check on the order e terms ef the left-hand-side of (23) which 
were originally obtained by means of tedious calculation. 
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The calculation of d' y/dt is straightforward. Second derivatives of F result 
but they can be eliminated by differentiating the lowest order form of (6) to find 


d' dF 


(MM + NN) : V(U + HL)— ^r- 

d ^ 


dF 


d * 


(f) 


- V • L 


dF 

dH 


(24a) 


^ H = -(£■ V)F + T V • £-7^- UL : V£fjf 


(24b) 


The final result for d' y/dt is 


e ^ 

e dt 


Lx = — f- 

it b|: 


• af : 7V 4 £ 5 * • If 


(t ♦ u • ar) ( 2aS • ar + 7 l t ■■ v£ ) 


,-dL 2 2 d , 

+ x’m- a^v • Lt-j-ar = vv 


l 2 

+ Th : VVI 


I 2 /„ dL .. 4 . _ \1 dF 

1 : VV_ T ( M ' dt + U<l : VL )| ~h 2 \ 

^ H*T/ 




’ dt 


A + UM:Vl)i a :W+I^ (m • fr 


<*» /»■ 


+ UH : VL 


) 
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Inserting the appropriate expressions (6) and (10) into the left-hand-side of (23) 
yields (25); which was to be shown. 

Given the results of (IQ) one may derive (I) for given (8) and constraint (4b), 

J ana can be obtained and thus (6). Some ambiguity is involved in the converse 
due to the arbitrariness involved in the choice of intermediate variables used to 
solve the single particle equation. 15 Nevertheless, a set of equations that is at 
least equivalent to those of (in) may be obtained by using (5) and (4b) to calcu- 
late transformation equations of the form of (8). These relations will yield J and 
y (though possibly not the same J and y) and (23) can be used to complete the 
derivation. The ambiguity is no shortcoming: those equations obtained could 
have been originally derived through the method of (HI). We simply remark that 
in using (HI) as a starting point the ambiguity is removed by specifying the al- 
lowed choice in variables. 
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Simplification in the form of each of these kinetic equations can be achieved 
by transforming velocity space variables so that the magnetic moment series 
makes up one component. The kinetic equation of (I) becomes 

dl 

Jt (r. M- V) + a(r, M. V) ’ Vf(r, /±, 77) + T(r, fi, 77) (r, /*, 77) = 0 (26) 

where 





(NM-MN) : V(U + tjL) 


d' fj. 

dt 




Note that f remains as originally defined, the probability density of particles 
in the r, <t 2 / 2, 77, v space, not in r, V, v ♦ 

In (IQ), anew distribution function f is defined such that ? represents the 
probability density of guiding centers in the space P, A, V| ( , $ where 


I 2 el „ d e 2* /3 ~ ~ 1 ~ 

A = IB" ' dt (U + HL)- g 2B : V ( U + HL -» ( 27 > 


and 
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The equation expressing conservation of guiding centers in this space is 


* & n - o. 


The egression for S ' P/S t can be obtained from (10a) by using (27) and 


= + U * aT + B AV • L 


+ eA \ l a tVLl^tVU- jl # : VUI^ : VL+I y : VUV-L 


d „ JL x (L 4 V) U 

+ LL: VUI y : VL + ^ l y : VL + BV- < § 


+ V • {NUM : VL-MUN: Vl} 


+ e V„ A I„ : VLV • L + BV 

ii y 


L x (L * V) L 


e dU ~ dL 
+ B l M * dt N ’ dt 


„ dU ~ dL\ 
- N • dfM * a T J 


In exactly the same manner as outlined above, we form the equations 


f d 3 PdAdV„ d4> = f d : 


'-(t) 


d r)d v 
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Gd 3 PdAdV N d$ = d 3 rd^— JdTjd 


and conclude that 


f(P. A, V,,) = Gf(r, fj., 7], v) = G 


|^f(P, A, 


- e A(NM-MN) : VL = G[? + y'] 


G = B 


f *1 


(NM-MK) : V(U + V H L) 


Inserting (29) into (27) one may obtain 


f ey 
G + G 




{« £ 


(NM-MN) : VL- A(NM-MN) : VLUL : VL> 


+ e A(L • V) f (31) 
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where 


D' _ d /S'P \ 8<v h d 

Dt at + \ St" ‘ y / + st av,. 


Substituting the quantities specified into (3 1) and performing the operations indicated 
yields equality. 

Since the sum of V • (8 ' P/St) and d/dV n (s ' V„ /8t)does not vanish. Equation 
(28) expresses the fact that D' f/Dt 4 0 and, like F, ? is not conserved along a 
guiding center trajectory. We postulate that a function that does satisfy such a 
requirement is K(P, A, V,, ) = f/G. We then form the expression: 


1 D* f __f_ D^_G 
G Dt “ G 2 Dt 


~± f / 8 *P\ d f S ' V ll\ 1 D* ol 
G [ V \St Y av (| \ St / + G Dt J 


= 0 . 


Thus, correct to first order in e , the function 


■ft {> + f 


(N M- MN) : V(U + 


v « £ >} 


is conserved along a guiding center trajectory. 

The additional results of (II) that are to be verified here are the first order 
corrections to the heat flux tensor Q. To follow the procedure of (DI) which was 
used to verify the pressure tensor of (II) (see footnote 5), we begin with the defi- 
nition of the heat flux tensor: 


fd 3 v(v-u)(v-u)(v-u) 
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Si i Slftw. iliiUfii li I 


where f is the solution of (1) and the definitions 


u 



(33a) 


and 


n 



(33b) 


will be used. 

In order to calculate Q correct to order e from the results of Section m, 
an equation of the general type (30a) must be enforced. Here we desire a differ- 
ent result from that obtained between f and f. This is the relation: 

fd 3 v = f J' d Ad V |( d $ (34a) 


where 


1 + e( 2 BA) 1/2 cos 


£(2 BA) 1/2 sin • (34b) 


Now using (34) and definitions (33) in (32) with v given by (9a) and (8) we can 
arrive at the following portions of Q previously presented in (II) 19 


r _ <s_ “ A 

[Q • L L] i ~ "g L x 


ps, --j-v-r* (R, - 3 Rj) (£ • V) L 
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+ 2 (Q„ “ Q) (L- V)« a + 2Q7u, | 

+ 2 q|mLM : Vu A + N L V : Vajfj 

\ [Q: (Ml8 + NN)] i = if Lx jvR a 


4P *1 

* ^(QtR^-R,} (£-V)£--/v-|>J 


where, to lowest order 


p = on 


u n = ij 


f V„ dAdV; 


r, = »ff(v„ - u„)>dA dv, ■ 

J 

II 

9 

j* f B A d A d V,j 

<3, = ”|f(vr u ii) SdAdV . • 

Q = m 

• 

ffEA(V,-u,) d A d Vj| 

R, = m r?(V,-u,)*dAdV, , 

Pd 

M 

II 

9 

« 

• 

f B A(V„ -u,,) 2 d A dV, 


R 3 = ml f [B A] 2 d Ad v ( | , 


P = P„ LL + P A (MM + Nl*) • 
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We hasten to note that one important result of (Q) that cannot be obtained 
from (HI) is that of the computation of u A correct to order e a . In principle, 
however, this can be accomplished by calculating the order e 2 portions of relations 
(8) and usin^ a higher order version of the procedure outlined above. 

To sum up, one must conclude that the computations provided here lend 
strong credibility to the assertion that all three of the papers treated are correct 
and equivalent. To be sure, one would expect this to be true, but given the vast 
observable differences pointed out in the introduction, and given the possibility of 
algebraic as well as theoretical error, this demonstration is altogether reassuring. 
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FOOTNOTES 


1. E. Frieman, R. Davidson, and B. Langdon, Phys. Fluids 9, 8 , 1475 (1966). 

2. A Macmahon, Phys. Fluids 8 , 10, 1840 (1966). 

3. G. Wilson, (to be published). 

4. This necessity is verified in (I) where a closed system of hydromagnetics 
equations is obtained in the two dimensional case only. 

5. See Ref. 3, eq. 38. 

6 . G. Chew, M, Goidberger, and F. Low, Proc. Roy. Soc. (London) A236 . 112 
(1956) 

7. This equation and 6 b are written in a slightly different form from that shown 
in (I) . For the purposes of comparison it is convenient to adhere strictly to 
the notion that f represents probability density in the velocity space 
element o-dod^dv and is, therefore, a function of a 2 /2 , 77 and v. 

■f 

8 . Ref. 1, pp. 1477-1479. I 

9. H. Grad, Proceedings of the Symposium on Electrodynamics and Fluid 
Dynamics of Gaseous Plasmas , (Polytechnic Press, Brooklyn, New York, 

1961, pp. 37-64). 

10. H. Grad, AEC Report NYO-1480-50, MF-48, New York University, 1966. 

11. Ref. 1, Sec. Q. 

12. M. Kruskal, J. Math, Phys. 3, 806 (1962). 

13. In principle. Relations ( 8 ) may be calculated to all orders in e. 

14. Ref. 11, p. 817, Sec. C.l. 
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15. See T. Northrop, M. Kruskal, and C. Liu, Phys. Fluids 9, 1503 (1966), 
pp. 1505-1506 for a discussion. 

16. R. Davidson, Phys. Fluids 10, 669 (1967). 

17. In (29) d/dt L Is the convective derivative of £, dL/dt + (S' P/8t * V)L cor- 
rect to order e. 

18. The fact that f/ B is conserved along the lowest order guiding center tra- 
jectory has been pointed out by Y. Whang (to be publicaed). 

19. Ref. 2, eqs. 14-16. Note a difference in definition by a factor of 2 between 
35a and that given in (II). 
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